Method of estimating initial values of potential in semiconductor device simulation

ABSTRACT

In order to estimate an initial potential value for semiconductor device simulation at each of iterative procedures a computer system, a plurality of bias conditions are stored in a memory. Following this, one bias condition is retrieved from the memory at a given iterative procedure. Further, an analysis result already obtained in an iterative procedure, which precedes the given iterative procedure, is retrieved from the memory. Subsequently, an initial potential value is estimated which is used in the give iterative procedure by solving a Laplace equation which is weighted by a coefficient including a reciprocal of electric field intensity.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The present invention relates generally to a computer device simulationand more specifically to a method of estimating initial values ofpotential during semiconductor device simulation.

2. Description of the Related Art

In the art of numerically analyzing a semiconductor device using acomputer, it is known to utilize drift-diffusion model which has beenapproximated by treating carriers (viz., electrons and holes) as afluid. Further, it is known in the art to use an energy transport modelestablished with higher order approximation compared with the abovementioned drift-diffusion model. It is to be noted that the presentinvention is applicable to both the aforesaid drift-diffusion model andthe energy transport model.

In a semiconductor device simulation using a drift-diffusion model undera stationary state, the following three fundamental equations (1), (5)and (6) should be established:

    divD=ρ                                                 (1)

    D=εE                                               (2)

    E=-gradψ                                               (3)

    ρ=q(p-n+N.sub.D -N.sub.A)                              (4)

where

D: electric flux density;

ρ: electric charge density;

E: electric field strength;

ε: dielectric constant;

ψ: potential;

q: elementary charge;

p: hole density;

n: electron density;

N_(D) : donor density; and

N_(A) : acceptor density.

    divJn=q·(R-G)                                     (5)

    divJp=-q·(R-G)                                    (6)

where

Jn: electron current;

Jp: hole current;

R: carrier recombination term; and

G: carrier generation term.

    Jn=q·n·μ.sub.n ·E+q·Dn·gradn(7)

    Jp=q·p·μ.sub.p ·E-q·Dp·gradp(8)

where

μ_(n) : electron mobility;

μ_(p) : hole mobility;

Dn: electron diffusion coefficient; and

Dp: hole diffusion coefficient.

    Dn=μ.sub.n ·{(k.sub.B ·T)/q}          (9)

    Dp=μ.sub.p ·{(k.sub.B ·T)/q}          (10)

where

k_(B) : Boltzmann's constant; and

T: temperature.

Equation (1) is called Poisson's equation. On the other hand, equation(5) is known as an electron current continuity equation while equation(6) is known as a hole current continuity equation.

In the above equations, the variable to be solved are the potential ψ,the electron density n, and the hole density p. Generally, a pluralityof input biases are firstly considered as boundary conditions. Then, thebasic equations (1), (5) and (6) are solved by sequentially renewing theinput biases.

Since each of the basic equations (1), (5) and (6) is a non-linearequation, it is solved using iterative calculations known as a Newtonmethod.

The above mentioned device simulation and the Newton method are known inthe art. By way of example, for further details reference should be madeto a book entitled "Very High Speed MOS Devices" edited by SusumuKohyama, published by Oxford University Press, 1990, pages 277-290.

In a device simulation, a two-dimensional region to be analyzed ispartitioned into a plurality of meshes. Subsequently, the abovementioned three basic equations (1), (5) and (6) are prepared orestablished at each mesh point. Therefore, assuming that N representsthe number of meshes, N variables develop with respect to each ofpotential, electron density, and hole density. Therefore, it isnecessary to solve 3N simultaneous equations.

Poisson's equation (1) is rewritten into

    F.sub.ψ (ψ,n,p)=0                                  (11)

The electron current continuity equation (5) is rewritten into

    F.sub.n (ψ,n,p)=0                                      (12)

Further, the hole current continuity equation (6) is rewritten into

    F.sub.p =(ψ,n,p)=0                                     (13)

In the above equations (11)-(13), the notations ψ, n, and p respectivelydenote potential, electron density, and hole density. Further, each ofthe notations ψ, n, and p represents N variables over an entiretwo-dimensional region. As is known in the art, two methods are proposedto solve the three basic equations (11)-(13). One method is asimultaneous (or coupled) method while the other method is Gummel'smethod or decoupled method. Although the simultaneous method is capableof solving the equations through a smaller number of iterativecalculations, it is sometimes unable to converge in the case where agood initial value is not given. On the contrary, the Gummel's methoddoes not rely on an initial value compared with the simultaneous method.However, according to the Gummel's method, a relatively large number ofiterative operations are required compared with the simultaneous method.As a result, in the event that a good initial value is given, thesimultaneous (i.e., couple) method is preferable over the Gummel'smethod in view of the fact that an overall calculation time can beshortened.

One known approach for obtaining initial potential values forimplementing iterative calculations of a device simulation, is proposedin a paper entitled "Initial guess strategy and linear algebratechniques for a coupled two-dimensional semiconductor equation solver"by Edwards, et al., NASECODE IV, 1985, pages 272-280.

However, the prior art proposed by Edwards, et al. has encountered thedifficulties that a large number of calculations are inevitablyrequired.

Therefore, what is required is a method of estimating initial values ofpotential with a small number of calculations.

SUMMARY OF THE INVENTION

It is the object of the present invention to provide a method ofestimating initial potential values during semiconductor devicesimulation with a small number of calculations.

An aspect of the present invention resides in a method of estimating aninitial potential value for semiconductor device simulation at each ofiterative procedures, said method being implemented using a computersystem, said method comprising the steps of: (a) storing a plurality ofbias conditions in a memory of said computer system; (b) retrieving onebias condition from said memory at a given iterative procedure; (c)retrieving an analysis result already obtained in an iterative procedurewhich precedes said given iterative procedure; and (d) estimating aninitial potential value for use in said given iterative procedure bysolving a Laplace equation which is weighted by a coefficient ω whichincludes a reciprocal of electric field intensity E.

BRIEF DESCRIPTION OF THE DRAWINGS

The features and advantages of the present invention will become moreclearly appreciated from the following description taken in conjunctionwith the accompanying drawings in which like elements are denoted bylike reference numerals and in which:

FIG. 1 is a flow chart which includes steps which characterize thepresent invention;

FIG. 2 is a cross section of an N-type MOSFET for use in discussing thepresent invention;

FIG. 3 is a graph showing an impurity profile at a portion of the deviceshown in FIG. 2;

FIG. 4 is a diagram showing a plurality triangular meshes formed on thecross section shown in FIG. 2;

FIG. 5 is a graph showing estimated potential values obtained using thepresent invention and a known method;

FIG. 6 is a cross section of an NPN-type bipolar transistor for furtherdiscussing the present invention;

FIG. 7 is a graph showing an impurity profile at a portion of the deviceshown in FIG. 6;

FIG. 8 is a diagram showing a plurality triangular meshes formed on thecross section shown in FIG. 6; and

FIG. 9 is a graph showing estimated potential values obtained using thepresent invention and a known method;

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

Before describing preferred embodiments of the present, it is deemedadvantageous to briefly discuss the principle underlying the presentinvention.

Let us consider that a series of bias conditions are sequentiallyanalyzed on a one-by-one basis. It is assumed that a given (e.g. ani-th) bias condition has been analyzed and the next (viz., (i+1)-th)bias condition is now to be analyzed. In such a case, Poisson's equationhas already been established in connection with the previous biascondition.

    divD.sub.i =ρ.sub.i                                    (14)

    D.sub.i =εE.sub.i =-εgradψ.sub.i       (15)

where

D: electric flux density;

ρ: electric charge density;

E: electric field intensity; and

ε: dielectric constant.

A bias potential of the current procedure is designated by ψ_(i+1)=ψ+δψ_(i) (where δψ_(i) represents the previous potential variation). Ifcarrier's re-distribution is not considered (viz., assuming that chargedensity's variation is not considered (i.e., ρ_(i+1) =ρ_(i))), then weobtain

    divD.sub.i =-εΔψ.sub.i =ρ.sub.i      (16)

    divD.sub.i+1 =-εΔ(ψ.sub.i +δψ.sub.i)=ρ.sub.i+1                        (17)

In other words, the above assumption implies that Poisson's equation(14) is only considered with both the electron and hole currentcontinuity equations neglected.

By subtracting the both sides of equation (16) from those of equation(17), the following Laplace equation is obtained

    εΔδψ.sub.i =0                      (18)

The value δψ which can be obtained by solving Laplace equation (18) isan approximate value of a potential variation. However, such a value isroughly approximated because re-distribution of carriers (viz.,electrons and holes) is not taken into account.

As is known, within a region in the vicinity of a PN junction, thecarrier density decreases due to recombination of electron-hole pairs.Further, a potential difference is large due to impurity ion charges insuch a region whereby high electric field intensity is generated. Thismeans that the resistance in such a region becomes high relative to theneighboring regions. Accordingly, the potential variations resultingfrom the bias applied is liable to propagate into the PN junctionregion. Therefore, in connection with Laplace equation (18), it can betreated that the dielectric constant is apparently small in the regionwhere high electric field intensity is exhibited.

Accordingly, if a coefficient of Laplace equation (18) is rendered smallin the region with high electric field intensity, it is possible topropagate more amount of potential variations into the high electricfield region. Thus, it can be regarded that the carrier re-distributionis desirably considered when solving Laplace equation (18). In otherwords, if a coefficient of Laplace equation (18) is a reciprocal of theelectric field intensity, a better approximate of the potentialvariation can be obtained. Therefore, it is possible to obtain betterinitial values of potential.

In such a case, the number of Laplace equations to be solved is equal tothe number of initial potential values to be determined and accordinglycan be halved compared with the conventional method.

In accordance with the above mentioned principle underlying the presentinvention, Laplace equation is solved which is weighted by a coefficientincluding a reciprocal of the electric field intensity as shown below.

    Δ(ω·δψ)=0                   (19)

    ω=(C.sub.s0 /C.sub.s1)exp(-α)                  (20)

    α={1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)}·1n{(E+C.sub.0)/C.sub.f1}             (21)

    ω=(C.sub.s0 /C.sub.s1){C.sub.f1 /(E+C.sub.0)}.sup.B  (22)

where

    B=1n(Cs0/Cs1)/1n(Cf0/Cf1)ψ.sub.i+1 =ψ.sub.i +δψ(23)

where:

δψ: potential variation;

ψ_(i) : potential which has been obtained by analyzing a previous biascondition;

ψ_(i+1) : initial value of potential which is used for a current biascondition;

E: electric field intensity;

C_(s0) : an empirically determined constant value corresponding to amaximum length characterizing a device configuration;

C_(s1) : an empirically determined constant value corresponding to aminimum length characterizing the device configuration;

C_(f0) : an empirically determined constant value corresponding to amaximum electric field intensity;

C_(f1) : an empirically determined constant value corresponding to aminimum electric field intensity; and

C₀ : a positive constant.

In the above, C₀ indicates a positive constant which ensure acalculation when the electric field intensity E is zero. Therefore, thepositive constant C₀ is appropriately selected considering the constantvalue C_(f1). On the other hand, the constant value C_(f1) is tonormalize the electric field intensity E and regulates the order of thecoefficient ω. The ratio of C_(s0) to C_(s1) (viz., C_(s0) /C_(s1))provides the maximum value in terms of the device's configuration to beanalyzed and also regulates the order of the coefficient ω. It isunderstood that the ratio of C_(f0) to C_(f1) (viz., C_(f0) /C_(f1))indicates the maximum variation of the electric field intensity.Finally, the logarithmic ratio {1n(C_(s0) /C_(s1))/1n(C_(f0) /C_(f1))}regulates coupling intensity between the variations of the electricfield and the coefficient ω.

FIG. 1 is a flow chart which depicts steps which characterize thepresent invention.

Firstly, a plurality of bias conditions to be analyzed are stored in anappropriate memory in a computer system (step 10). By way of example,the number of bias conditions is denoted by N. Subsequently, at step 12,a counter is reset which is provided for counting the number of biasconditions which have been analyzed. At step 14, the content of thecounter is incremented by one. At the first iterative procedure, thefirst bias condition is retrieved from the memory at step 16, afterwhich the routine goes to step 18 at which a check is made to determineif the content of the counter equals one (1). In this case (viz., thefirst iteration), the answer is affirmative and thus, the routine jumpsto step 24 at which the above mentioned three device equations (1), (5)and (6) are solved using a simultaneous (or coupled) method. Followingthis, the result obtained by the analysis implemented at step 24 isstored in the memory (step 26), and the program goes to step 28. In thiscase, the content of the counter does not reach N and thus the routinereturns to step 14.

If the outcome of the inquiry performed in step 18 is "NO", the programproceeds to step 20 at which the previously obtained analysis result isretrieved from the memory. Following this, at step 22, an initial valueof potential is estimated by solving the Laplace equation (19) which isweighted by a coefficient which is a reciprocal of the electric fieldintensity (E) as above mentioned. Thereafter, the routine goes to step24.

In the event that the answer at step 28 is "YES", the programterminates.

The present invention is further discussed in detail with two computersimulations which were conducted by the inventor. Firstly, the presentinvention was applied to an N-type MOSFET (Metal Oxide-SemiconductorField-Effect Transistor), which is discussed with reference to FIGS.2-5.

FIG. 2 is a diagram schematically showing a cross section of an N-typeMOSFET. The configuration of the device shown in FIG. 2 is readilyunderstood and hence no further description is given for the sake ofsimplifying the disclosure.

FIG. 3 is a graph showing an impurity profile of the device of FIG. 2 ata depth of Y=1.0 μm (see FIG. 2).

FIG. 4 is a diagram showing the manner wherein the cross section of FIG.2 is partitioned into 1475 triangular meshes for obtaining potentialvalues shown in FIG. 5.

FIG. 5 is a graph showing the relationship between the two estimatedinitial values of potential (denoted by curves B and C) and thepotential obtained by carrying out using the just mentioned initialvalues (denoted by A). More specifically, the curve B indicates theestimated initial value of potential obtained using the presentinvention, while the curve C indicates the estimated initial value ofpotential obtained using the above mentioned prior art (viz., the methodproposed by S. P. Edwards, et al).

In order to obtain the relationship shown in FIG. 5, the bias conditionat the previous iterative calculation was set such that 1.0V was appliedto the gate electrode, 0.5V being applied to the drain electrode, and0.0V being applied to both the source and substrate electrodes (FIG. 2).At the next (viz., current) iteration, only the drain bias was increasedby 0.1V. It is understood that the estimated initial value of potentialaccording to the present invention is nearer to the analysis resultindicated by the curve A compared with the estimated initial potentialobtained by the known technique.

The number of iterative calculations for estimating the initialpotential of each mesh was four (4) according to the present inventionwhilst amounting to five (5) according to the known technique. In thiscase, the converging condition was set to 1.0×10⁻⁶. Further, with acomputer equipped with a CPU (Central Processing Unit) of 33 MIPS(Million Instructions Per Second), the time required for estimating eachof the above mentioned initial potential values was 4.35×10⁻³ secondsaccording to the present invention. On the contrary, according to theprior art, 6.23×10⁻³ seconds was necessary for the same goal.

Secondly, the present invention was applied to an NPN-type bipolartransistor, which is discussed with reference to FIGS. 6-9.

FIG. 6 is a diagram schematically showing a cross section of theNPN-type bipolar transistor. The configuration of the device shown inFIG. 6 is readily understood and hence no further description is notgiven for the sake of simplifying the disclosure.

FIG. 7 is a graph showing an impurity profile of the device of FIG. 6 ata location of X=2.55 μm (see FIG. 6).

FIG. 8 is a diagram showing the manner wherein the cross section of FIG.6 is partitioned into 1495 triangular meshes for obtain potential valuesshown in FIG. 9.

FIG. 9 is a graph showing the relationship between the two estimatedinitial values of potential (denoted by curves B and C) and thepotential obtained by carrying out using the just mentioned initialvalues (denoted by A). More specifically, the curve B indicates theestimated initial value of potential obtained using the presentinvention, while the curve C indicates the estimated initial value ofpotential obtained using the above mentioned prior art (viz., the methodproposed by S. P. Edwards, et al).

In order to obtain the relationship shown in FIG. 9, the bias conditionat the previous iterative calculation was set such that 0.8V was appliedto the base electrode, 1.05V being applied to the collector electrode,and 0.0V being applied to both the emitter and substrate electrodes (seeFIG. 6). At the next (viz., current) iteration, only the base bias wasincreased by 0.02V. It is understood that the estimated initial value ofpotential according to the present invention is nearer to the analysisresult indicated by the curve A compared with the estimated initialpotential obtained by the known technique.

The number of iterative calculations for estimating the initialpotential of each mesh was five (4) according to the present inventionwhilst amounting to six (5) according to the known technique. In thiscase, the converging condition was set to 1.0×10⁻¹². Further, with acomputer equipped with a CPU (Central Processing Unit) of 33 MIPS(million Instruction Per Second), the time required for estimating eachof the above mentioned initial potential values was 5.26×10⁻³ secondsaccording to the present invention. On the contrary, according to theprior art, 5.29×10⁻³ seconds was necessary for the same goal.

From the foregoing, it is understood that the present invention is ableto provide good initial values of potential without rendering complexthe numerical calculation. Therefore, the present invention is able toeffectively shorten the calculation time for estimating initialpotential values for implementing potential analysis in a devicesimulation.

It will be understood that the above disclosure is representative ofonly few possible embodiments of the present invention and that theconcept on which the invention is based is not specifically limitedthereto.

What is claimed is:
 1. A computer-implemented semiconductor devicesimulation method comprising the steps of:(a) partitioning atwo-dimensional analysis region of the semiconductor device into Nmeshes and storing N bias conditions respectively associated with said Nmeshes; (b) generating numerical values used in the semiconductor devicesimulation by solving a set of equations representative of saidsemiconductor device by the following substeps including:(1) retrievinga first of said N bias conditions; (2) solving said set of equationsusing said first bias condition; and thereafter (3) solving said set ofequations using one of remaining N-1 bias conditions and an initialpotential value ψ_(i+1) which is estimated by solving a Laplace equationdefined by:

    Δ(ω·δψ)=0

where

    ω=(C.sub.s0 /C.sub.s1)exp(-α)

    α={1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)}·1n{(E+C.sub.0)/C.sub.f1 }

such that

    ω=(C.sub.s0 /C.sub.s1){C.sub.f1 /(E+C.sub.0)}.sup.B,

where

    B=1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)ψ.sub.i=1 =ψ.sub.i +δψ

and where δψ: potential variation; ψ_(i) : potential obtained from adirectly-preceding solution of said set of equations; E: electric fieldintensity; C_(s0) : a constant value corresponding to a maximum lengthcharacterizing a device configuration; C_(s1) : a constant valuecorresponding to a minimum length characterizing the deviceconfiguration; C_(f0) : a constant value corresponding to a maximumelectric field intensity; C_(f1) : a constant value corresponding to aminimum electric field intensity; and C₀ : a positive constant; and(4)consecutively repeating substep (3) for remaining N-2 bias conditions;and (c) simulating the semiconductor device using said numerical values.2. The method according to claim 1, wherein the set of equationscorresponds to a Poisson's equation.
 3. The method according to claim 1,wherein the set of equations is solved using a simultaneous method. 4.The method according to claim 1, wherein the set of equations is solvedusing a coupled method.
 5. A computer-implemented semiconductor devicesimulation method comprising the steps of:(a) partitioning an analysisregion of the semiconductor device into N analyzing regions and storingN bias conditions respectively associated with said N analyzing regions;(b) retrieving a first of said N bias conditions; (c) generating andstoring a numerical analysis result for an analyzing region associatedwith said first bias condition based on said first bias condition; (d)generating and storing numerical analysis results for remaining N-1analyzing regions on a one-by-one basis based on remaining N-1 biasconditions, wherein an initial potential value ψ_(i+1) is used ingenerating the numerical analysis results for said remaining N-1analyzing regions, said initial potential value ψ_(i+1) being estimatedfor each of said remaining N-1 analyzing regions by solving a Laplaceequation defined by:

    Δ(ω·δψ)=0

where

    ω=(C.sub.s0 /C.sub.s1)exp(-α)

    α={1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)}·1n{(E+C.sub.0)/C.sub.f1 }

such that

    ω=(C.sub.s0 /C.sub.s1){C.sub.f1 /(E+C.sub.0)}.sup.B,

where

    B=1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)ψ.sub.i=1 =ψ.sub.i +δψ

and where δψ: potential variation; ψ_(i) : potential obtained from adirectly-preceding numerical analysis result; E: electric fieldintensity; C_(s0) : a constant value corresponding to a maximum lengthcharacterizing a device configuration; C_(s1) : a constant valuecorresponding to a minimum length characterizing the deviceconfiguration; C_(f0) : a constant value corresponding to a maximumelectric field intensity; C_(f1) : a constant value corresponding to aminimum electric field intensity; and C₀ : a positive constant; and (e)simulating the semiconductor device using said numerical values.
 6. Themethod according to claim 5, wherein said N bias conditions represent aset of boundary conditions at each of said N analyzing regions.
 7. Themethod according to claim 6, wherein the step of generating a numericalanalysis result in steps (c) and (d) includes the step of solving aPoisson's equation.
 8. The method according to claim 6, wherein the stepof generating a numerical analysis result in steps (c) and (d) includesthe step of solving a set of device equations using a simultaneousmethod.
 9. The method according to claim 6, wherein the step ofgenerating a numerical analysis result in steps (c) and (d) includes thestep of solving a set of device equations using a coupled method. 10.The method according to claim 5, wherein said N analyzing regionscomprise of N triangular meshes.
 11. A computer used in semiconductordevice simulation, said computer being programmed to carry out the stepsof:(a) partitioning a two-dimensional analysis region of thesemiconductor device into N meshes and storing N bias conditionsrespectively associated with said N meshes; (b) generating numericalvalues used in the semiconductor device simulation by solving a set ofequations representative of said semiconductor device by the followingsubsteps including:(1) retrieving a first of said N bias conditions; (2)solving said set of equations using said first bias condition; andthereafter (3) solving said set of equations using one of remaining N-1bias conditions and an initial potential value ψ_(i+1) which isestimated by solving a Laplace equation defined by:

    Δ(ω·δψ)=0

where

    ω=(C.sub.s0 /C.sub.s1)exp(-α)

    α={1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)}·1n{(E+C.sub.0)/C.sub.f1 }

such that

    ω=(C.sub.s0 /C.sub.s1){C.sub.f1 /(E+C.sub.0)}.sup.B,

where

    B=1n(C.sub.s0 /C.sub.s1)/1n(C.sub.f0 /C.sub.f1)ψ.sub.i=1 =ψ.sub.i +δ104

and where δψ: potential variation; ψ_(i) : potential obtained from adirectly-preceding solution of said set of equations; E: electric fieldintensity; C_(s0) : a constant value corresponding to a maximum lengthcharacterizing a device configuration; C_(s1) : a constant valuecorresponding to a minimum length characterizing the deviceconfiguration; C_(f0) : a constant value corresponding to a maximumelectric field intensity; C_(f1) : a constant value corresponding to aminimum electric field intensity; and C₀ : a positive constant; and(4)consecutively repeating substep (3) for remaining N-2 bias conditions;and (c) simulating the semiconductor device using said numerical values.12. The computer according to claim 11, wherein said computer has amemory, and said N bias conditions and solutions to said set ofequations are stored in said memory.
 13. The computer according to claim12, wherein said N bias conditions represent a set of boundaryconditions at each of said N meshes.